#include "bc_constantWaterDepth_velocity.hpp"
#include "meshTypes.hpp"
#include "field.hpp"
#include "boundaryPatch.hpp"
#include "sweFunctions.hpp"
#include "simulationTime.hpp"
#include "exceptions.hpp"
#include "geometry.hpp"
#include "sweSolver.hpp"

namespace hms {

#define CHEAD(RETURN_TYPE) \
template<typename MeshType> \
RETURN_TYPE ConstantWaterDepth_velocity<MeshType>::

CHEAD() ConstantWaterDepth_velocity(
	Field<MeshType>& v,
	const BoundaryPatch& patch,
	const Field<MeshType>& dFieldArg,
	scalar dGhostArg
) :
	BaseType( Type::constantWaterDepth_velocity, v, patch ),
	m_transformer { v, patch }
{
	dField( dFieldArg );
	dGhost( dGhostArg );
}

CHEAD(void) update() const {
	Field<MeshType>& v { field() };
	const Field<MeshType>& d { dField() };
	/* eigen conditional type. Basically cell-wise true/false with true if
	 * supercritical or critical (Fr >= 1) */
	auto supercritical {
		isSupercritical( v.domainValues( patch() ), d.domainValues( patch() ) )
	};

	Array2Xs& vTransformed { m_transformer.transform() };

	/* supercritical INflow is a problem, outflow is not */
	if ( (supercritical && vTransformed.row(0) < 0).any() )
		throw BoundaryConditionError(
			"Supercritical inflow detected at boundary patch \""
			+ patch().name() + "\". "
			"Boundary condition \"ConstantWaterDepth\" cannot be used. "
			"For supercritical inflow conditions, all state variables have "
			"to be defined by value."
		);
	
	/* supercritical outflow -> vGhost_ortho = vDomain_ortho
	 * subcritical in-/outflow ->
	 * vGhost_ortho = vDomain_ortho + 2 * sqrt(g*d_Domain) - 2 * sqrt(g*d_ghost)
	 */
	vTransformed.row(0) = supercritical.select(
		vTransformed.row(0),
		vTransformed.row(0) + 2*(
			sqrt(phys::g * d.domainValues( patch() )) - 
			sqrt(phys::g * m_dGhost) 
		)
	);

	v.ghostValues( patch() ) = m_transformer.retransform();
}

CHEAD(scalar) dGhost() const    { return m_dGhost; }
CHEAD(void  ) dGhost(scalar arg){ m_dGhost = arg; }

#undef CHEAD

#define INSTANTIATE_FOR_ALL_MESHTYPES(MeshType) \
template class ConstantWaterDepth_velocity<MeshType>;
#include "meshInstantiator.hpp"


} // end namespace hms